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Abstract 

High parton densities in ultra-relativistic nuclear collisions suggest a description of these 
collisions wherein the high energy nuclear wavefunctions and the initial stages of the nuclear 
collision are dominated by classical fields. This underlying paradigm can be significantly im- 
proved by including quantum fluctuations around the classical background fields. One class 
of these contributes to the energy evolution of multi-parton correlators in the nuclear wave- 
functions. Another dominant class of unstable quantum fluctuations grow rapidly with proper 
time r after the collision. These secular terms appear at each loop order; the leading contri- 
butions can be resummed to all loop orders to obtain expressions for final state observables. 
The all-order result can be expressed in terms of the spectrum of fluctuations on the initial 
proper time surface. We compute, in A T — gauge, the essential elements in this fluctuation 
spectrum-the small quantum fluctuation modes in the classical background field. With our 
derivation in QCD, we have all the ingredients to compute inclusive quantities in heavy ion 
collisions at early times including i) all-order leading logs in Bjorken 2:1,2 of the two nuclei, 
ii) all strong multiple scattering contributions, and iii) all-order leading secular terms. In the 
simpler analogous formalism for a scalar 4> 4 theory, numerical analysis of the behavior of the 
energy-momentum tensor is strongly suggestive of early hydrodynamic flow in the system [1] . 
In QCD, in addition to studying the possible early onset of hydrodynamic behavior, additional 
important applications of our results include a) the computation of sphaleron transitions off- 
equilibrium, and b) "jet quenching" , or medium modification of parton spectra, in strong color 
fields at early times. 

1 Introduction 

The large flow measured in heavy ion collisions at RHIC [2-5] and more recently at the LHC [6] 
can be described in hydrodynamic models that have both a nearly perfect fluid value of the shear 
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viscosity to entropy ratio of the quark-gluon matter produced and fairly short thermalization times 
that usually range between 0.5 and 2 fermis/c [7-9] (depending on the assumptions made about the 
initial conditions and the implementation of the freeze-out). How isotropization and (subsequently) 
thermalization is achieved in heavy ion collisions is an outstanding problem which requires that one 
understand the properties of the relevant degrees of freedom in the nuclear wavefunctions and how 
these degrees of freedom are released in a collision to produce quark-gluon matter. An ab initio 
approach to the problem can be formulated within the Color Glass Condensate (CGC) effective 
field theory, which describes the relevant degrees of freedom in the nuclei as dynamical gauge 
fields coupled to static color sources [10,11]. The computational power of this effective theory is 
a consequence of the dynamical generation of semi-hard saturation scale [12,13] larger than the 
intrinsic non-perturbative QCD scale {Q 2 S » Aq CD ), which allows for a weak coupling treatment 
of the relevant degrees of freedom [14-16] in the high energy nuclear wavefunctions. 

There has been significant recent progress in applying the CGC effective field theory to studying 
the early time behavior of the quark-gluon matter called Glasma [17] produced in the initial little 
bang of a high energy heavy ion collision. Inclusive quantities such as the pressure and the energy 
density in the Glasma can be written as expressions that factorize, to leading logarithmic accuracy 
in the longitudinal momentum fraction x, the universal properties of the nuclear wavefunctions 
(measurable for instance in proton-nucleus or electron-nucleus collisions) from the final state evo- 
lution of the matter in the collision [18-20]. Key to this approach are the quantum fluctuations 
around the classical fields. In particular, quantum fluctuations that are invariant under boosts can 
be shown to factorize into universal density functionals that encode the multi-parton correlations 
in the nuclear wavefunctions. The evolution of these density functionals with energy is described 
by the JIMWLK renormalization group equation [21-28]. 

There are however quantum fluctuations that are not boost invariant. It was observed in [29- 
31], via numerical solutions (see also [32,33] for a semi-analytic discussion of some instabilities in 
the solutions Yang-Mills equations) of the classical Yang-Mills equations, that rapidity dependent 
quantum fluctuations in the expanding Glasma are unstable and grow exponentially as the square 
root of the proper time r after the collision. In fact, both the existence and the specific time 
dependence of these instabilities was anticipated based on studies of the Weibel instabilities in 
expanding anisotropic Yang-Mills plasmas [34-37]. The unstable quantum fluctuations (initially of 
order 0(1)) become comparable in size to the classical field (of order 0(g~ 1 )) on a very short time 
scale t ~ Qj 1 . Fortunately, one can isolate and resum these rapidly growing secular divergences 
to all orders in perturbation theory. The resulting expressions are free of secular divergences, and 
can be rephrased as an average over a spectrum of Gaussian fluctuations of the initial data for 
the classical field encountered at leading order. A similar observation was made previously in the 
context of inflationary cosmology [38,39]. 

In a previous paper [1], we developed this formalism for a scalar <p 4 theory which, like QCD, has 
a dimensionless coupling in 3+1 dimensions and has unstable modes. We computed the spectrum 
of fluctuations and showed that the resummed expression for the pressure and energy density can 
similarly be expressed as an ensemble average over quantum fluctuations. The rapid growth of 
the unstable fluctuations has drastic consequences. Without resummation, the relation between 
the energy density and the pressure is not single valued. For the resummed expressions, while the 
relation between the pressure and energy density is not single valued at early times, it becomes so 
after a finite time evolution. This development of an equation of state therefore allows one to write 
the conservation equation for the resummed energy momentum tensor T' 1 " as the closed form set of 
equations corresponding to the hydrodynamical evolution of a relativistic fluid. This result can be 
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interpreted as arising from a phase decohcrence of the classical field configurations with different 
initial conditions given by the ensemble of quantum fluctuations. In this theory, the period of the 
classical trajectories is proportional to the amplitude of the field. Anharmonicity occurs in any 
non-linear system and we expect the same to hold for QCD. As the different trajectories become 
phase shifted for different amplitudes there are cancellations resulting in a single valued relation 
between the pressure and the energy density. 

This phenomenon shares several common features with Srednicki's hypothesis of eigenstate 
thermalization and Berry's conjecture [40-44]. Berry conjectured in [40] that high lying energy 
eigenstates of systems whose classical counterpart is chaotic have very complicated wavefunctions 
that for many purposes behave like random Gaussian 1 functions. A system in such an eigen- 
state would display features reminiscent of thermal equilibrium, despite being in a pure quantum 
state [42]. For a system starting initially in a coherent state rather than an energy eigenstate, 
thermalization would merely amount to losing the initial coherence. Although these ideas where 
formulated in much simpler systems, they may have some relevance to QCD since here also the 
underlying classical theory is believed to be chaotic [45,46]. 

In this paper, we shall focus on computing the initial spectrum of fluctuations in the Glasma 
formed at early times after a heavy ion collision. The classical background field at r = + in 
the Glasma can be expressed, from the continuity of the Yang-Mills equations across the light- 
cone [47,48], in terms of classical solutions of the Yang-Mills equations for each of the two nuclei 
before the collision. For later times, analytical solutions are not known 2 ; however, the Yang-Mills 
equations have been solved numerically with the initial conditions at r = + [51-55]; for a nice 
review, see [56]. Fortunately, inclusive quantities such as components of the energy-momentum 
tensor are sensitive only to the initial spectrum of fluctuations about the classical field at r = + , 
which can be calculated analytically. Specifically, we will solve the small fluctuations equations of 
motion in A T = gauge, in order to obtain a complete orthonormal basis of these fluctuations. 
There was a first attempt to compute the small fluctuations 2-point correlator in the Glasma [57] 
which, as we shall discuss, was incomplete because it did not include fully the structure of the 
background field. 

The paper is organized as follows. In the next section, we will outline the power counting 
of higher order contributions in the Glasma and emphasize the necessity of resumming secular 
terms. We isolate the leading contributions and obtain an expression for the resummed leading 
secular divergences. We show that this expression for inclusive quantities can be rewritten as a 
path integral over a spectrum of fluctuations times the leading order (classical) expression for the 
inclusive quantity. The only unknown ingredient in this reformulation are the small fluctuation 
fields on the initial proper time hypersurface. Additional sub-sections discuss gauge invariance 
issues and the renormalization of ultraviolet divergences. In section 3, we will show how to compute 
the small fluctuation fields in the vacuum. We first obtain an inner product for fluctuations on 
a space-like Cauchy surface that we use to define the orthogonality between a pair of fields. We 
will further prove that the inner product is independent of the chosen surface. We then show that 
the small fluctuation fields can be expressed as a linear combination of modes whose coefficients 
are Gaussian-distributed random complex numbers. (This is equivalent to diagonalizing the small 
fluctuation correlator on the initial proper time surface.) Because even the computation of the 

1 Naturally, the wavefunction of a given eigenstate is not a random function. Berry's conjecture means that 
for the purpose of computing the expectation value of sufficiently inclusive observablcs, one can replace the true 
wavefunction by a random Gaussian function. 

2 For some interesting recent attempts, see [49,50]. 
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free fluctuations in A T = gauge is non-trivial, we shall first solve small fluctuation equations 
in the vacuum. Then, in the section 4, we shall solve the small fluctuation equations in the 
background classical field of the Glasma to construct the corresponding physical small fluctuation 
modes. Section 5 outlines a practical algorithm to compute inclusive quantities (including all 
leading logarithms in x, and all leading secular contributions) as a function of proper time. As 
we shall demonstrate, the complexity of this space-time evolution is manageable, and amounts to 
diagonalizing certain matrices on the initial proper time surface. In the final section, we re-state 
our key results and discuss some important applications. These include a) a systematic study 
of possible thermalization of the quantum system whose evolution with proper time we plan to 
simulate numerically just as for the scalar case studied previously, b) The nature and role of 
sphaleron transitions in the early time dynamics of the system, c) The medium modification of 
hard probes to study 'jet quenching' at early times. An open question for future work is to explore 
until what times these analysis are valid and how one can incorporate sub-leading contributions 
that become increasingly important at late times. There are two appendices. The first discusses 
technical aspects of the computation of small fluctuation fields in the vacuum. Expressions for the 
Wightman functions for free fields in A T = gauge are discussed in the second appendix, where 
some connections to previous work on these is also discussed [57,58]. 



2 Resummation of leading instabilities in the Glasma 

We will begin by first outlining how the power counting for computing inclusive quantities in field 
theories with strong time dependent sources is modified due to the presence of secular divergences. 
Following this power counting, we derive an explicit expression for the energy-momentum tensor in 
heavy ion collisions that resums the leading instabilities to all loop orders in perturbation theory. 
We will show that the resummed expression for the energy momentum tensor can be expressed 
as a path integral over the product of two terms. The first is a weight functional that samples 
the spectrum of quantum fluctuations on the initial proper time hypersurface, while the second is 
the leading order expression for the energy-momentum tensor. In the latter, the classical field is 
shifted by the sampled quantum fluctuations. Computing the initial spectrum of fluctuations is 
our primary goal in this paper, the derivation of which will be discussed at length in sections 3 
and 4. Before we go there, two further sub-sections will discuss the constraints imposed by gauge 
invariance on the spectrum of fluctuations and the nature of ultraviolet divergences respectively. 



2.1 Power counting of unstable modes in the Glasma 

In previous works [59,60] , it was shown that the problem of computing leading order (LO) and next- 
to-leading order (NLO) contributions to inclusive quantities -such as components of the energy 
momentum tensor- in field theories with strong time dependent sources can be formulated as an 
initial value problem where a classical field determined on an initial Cauchy surface is evolved up 
to the time at which the (local) observable is computed. Because one anticipates that a semi- 
hard scale Q 2 S 3> A^ CD is generated by the non-linear QCD dynamics at high energy [12,13], a 
systematic weak coupling expansion of these inclusive quantities is feasible. One can formally 
arrange the perturbative expansion of an observable such as the energy momentum tensor as 

0[p!,p 2 } = ^ co + ag 2 + c 2 g 4 -\ , (1) 
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where each term corresponds to a different loop order. Each of the coefficients c n is itself an infinite 
series of terms involving arbitrary orders in (gpi, 2 ) F '■ These terms are all of order unity because 
the color charge densities are of order pi >2 ~ 0(g^ 1 ) in a large nucleus at high energy. The color 
charge densities correspond to the large x color sources in either nucleus 1 or nucleus 2 respectively 
in a heavy ion collision. Their evolution with the separation scale between sources and fields is 
described by the JIMWLK equation, which will be stated shortly. The LO contribution comes 
from the first coefficient Co, 

O lo [ PuP2 }^^. (2) 

This leading term co/g 2 has been studied extensively for the single inclusive gluon distribution in 
A+A collisions [51-54] and recently for the double inclusive distribution as well [61]. 
Following this terminology, we denote 

°nlo[Pi>P2] = ci , N ^ [pi,p 2 ] = c 2 g 2 , ••• (3) 

At each order in the loop expansion, there can arise contributions from the loop integrals which 
are of the same magnitude as lower orders. One set of such contributions are the increasingly 
large logarithms of the momentum fractions 2:1,2 of partons in the nuclear wavefunctions as higher 
energies, or equivalently smaller values of 0:1,2, are achieved in nuclear collisions. The term c„ can 
have up to n powers of such logarithms, with leading logarithmic terms identified as terms that 
have as many logarithms as their order in the loop expansion, 



1 00 ( \ 

LL o g [pi,P2] = -2 E dn [a 2 ln ( — 

• y n=0 K 12 



(4) 



where d n is the coefficient of n-th term in the leading log expansion. We were able to show [18-20] 
that the leading logarithmic contributions in x li2 , after averaging over the sources pi, 2 factorize 
into the expression 

(0> LLog = J [D Pl D P2 ] W X1 [pi] W X2 [ P2 ] LO [ Pl , P2 ] , (5) 

where W Xl 2 [pi t2 ] are the density functional we alluded to previously. These obey the JIMWLK 
equation [21-28] 

dW Xl2 [p h2 ] 
91n(l/a;i,2) 

Here ^1,2 are the JIMWLK Hamiltonians of the two nuclei; since their explicit form is not essential 
to the discussion here, we will refer the interested reader to ref. [18] for explicit expressions in our 
notation. Given an initial condition at some initial Xo value, the JIMWLK equation describes the 
evolution in the nuclear wavefunctions of the multi-parton correlators that contribute to inclusive 
observables measured in the final state. 

The resummation of quantum corrections arising from logarithms in £1,2, as sketched here, 
takes into account contributions that are essential in describing the energy evolution of inclusive 
observables in heavy ion collisions. These contributions are zero modes in u, the Fourier conjugate 
of the space-time rapidity n, and are localized in rapidity around the wave functions of the incoming 
nuclei. There are also quantum fluctuations that are non-zero modes of v. Such contributions, that 
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do not bring leading logs of l/xi.2, cannot be factorized into the evolution of the density functionals 
W Xl 2 in eq. (5). As shown in [29-31,62], these terms can be unstable and grow exponentially with 
the square root of the proper time (equal to te V2x + x~ in light-cone co-ordinates) for a system 
undergoing one dimensional longitudinal expansion. Based on these considerations, the expansion 
we sketched in eqs. (1) and (4) needs to be modified to keep track also of quantum fluctuations 
of amplitude gexp(y/jrF) (where \i is a growth rate of order Q s ) relative to the leading term. 
This is necessary because the rapid growth of these unstable modes leads to a break down of the 
perturbative expansion when 

t - r max ee ^ In 2 (7) 

is reached, the proper time at which 1-loop corrections become as large as the leading order term. 
The breakdown of the expansion can be avoided if one rcsums these divergent contributions, leading 
to a resummed result that is well behaved for t — > +oo. Taking into account both the leading logs 
in l/xi.2 and the leading unstable contributions, the new expansion reads 



O LLog [ Pl ,p 2 ] ee 1 £ g 2n £ d p , q In* (—) e 2 «V^ . 



(8) 



Thus far, we have only resummed the q = sector of this formula, where the result of the resumma- 
tion is expressed by the factorized formula (5). The two sources of leading quantum fluctuations at 
this accuracy can be resummed independently because a given quantum fluctuation mode cannot 
be at the same time a zero mode (that generates logarithms in 2:1,2) and a non-zero mode (that 
generates a secular divergence in proper time r). Naturally, in higher loop corrections, one loop can 
bring a log of 1/2:1,2 while another loop brings a secular divergence. This is why eq. (8) has terms 
with both p and q non-zero simultaneously. But the independence of the two types of divergences, 
based ultimately on considerations of causality, leads us to expect that the double series of eq. (8) 
can be factorized into a series in p times a series in q. 

2.2 All orders resummation of the leading secular terms 

We shall now discuss how resumming the leading secular terms modifies the expression of eq. (5). 
Albeit our considerations apply to any inclusive quantity, for specificity, we shall consider here the 
energy-momentum tensor. 

2.2.1 Reminder of LO and NLO results 

Let us recall first that at leading order in g 2 , the energy- momentum tensor T££ is given by 

= \^ v K^a^ - JTFL* , (9) 

with the field strength tensor defined as 

JT = d»Al - d v A» + gf abc A^A v c , (10) 
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where A% is the solution of the classical Yang-Mills 3 equations with sources that vanishes at 

x° — > — oo, 

VfT» v = 5 v+ p a l +5»-p a 2 , lim A%(x)=0. (11) 

One can then express the NLO contribution to the energy-momentum tensor, for a given distri- 
bution of color sources 4 and at an arbitrary space-time point, as the action of a functional operator 
acting on the LO result [18,59], 

T £o( x ) = [/ d 3 u(3-T u + ^J d 3 u d 3 v J2 J jSjkk [a+kXa ' T «H«-fcAa • T -l] . (12) 

where E is a Cauchy surface where the initial values of the classical color field and its derivatives 
are specified. In this formula, ft and a±Afea & r e small corrections to the gauge field A^. More 
specifically, (i(u) is the one loop correction to the classical field on the surface E (sec [18] for 
more details). In applications to heavy ion collisions, a natural choice of £ is the surface at 
proper time r = + , which corresponds physically to times just after the two nuclei have collided. 
Though t = + is the initial surface of choice, we will at the outset consider a generic space-time 
hypersurface. The only constraint on £ is that, for the forthcoming resummation to be effective, it 
must be located at times before the unstable modes have become too large 5 . In the right hand side 
of eq. (12), T u is the generator of shifts of the initial data for the classical field on E. Generically, 
it reads 6 

°- T »="^i^r {dViun miwr <13) 

where A^(u) (in curly font and without a time argument) is the value of the classical field on 
the initial time surface. Note that in general, specific gauge conditions and specific choices of the 
surface E reduce the number of terms this operator contains. In its minimal form, it contains one 
term for each independent field component or field derivative component that one must specify in 
the initial value problem on E. 

Let us now explicit a bit more the fields a^fcAa that appear in cq. (12). They are small 
fluctuation fields about the classical field A 11 , that obey the equation of motion 7 

£V (W*a v - V v a») - igF u >" = , (14) 



3 We have written the Yang-Mills equations in a form that involves the adjoint representation of the covariant 
derivative, T>® b = c9 M <5 ab — igA® b , where A ab is the classical gauge potential in the adjoint representation. It 
is important to distinguish the Aff's from the .4^'s that are the components of the SU(3) element _4 M in its 
decomposition over the generators of the algebra, _4 M = A^t a . The two sets of coefficients are related by A b ^ = 
~if abc Af l , since the components of the generators in the adjoint representation are (*2dj ) bc = ~if abc - 

4 Unless specified otherwise, the dependence on pi,p2, the color charge densities in each of the nuclei, will be 
implicit in our discussion. 

5 The NLO expression in eq. (12) does not depend on the choice of E. However, our resummed result will depend 
on this choice since it includes only a subset of the higher loop corrections. Provided the surface £ is located 
in a region where the unstable fluctuations are still small, the difference between various choices of £ is a small 
correction. We will discuss this point further later in the paper. 

6 It has dimension of (mass) 2 because dim[^] = (mass) 2 and dim[a-|-fcAa] = (mass) . 

7 To avoid cumbersome notations, we have not written explicitly the color indices of the various objects. Here, 
and henceforth, the T>'s and T should be understood as objects in the adjoint representation. For instance T v ' M a M 
Likewise, V ti V^a v is V ahzy .V^.a v c . 
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and that are specified in the remote past by the boundary condition 

lim a^ a (x)=^r e ^. (15) 

The T a 's are the SU(3) generators and e^ x is the polarization vector. Thus the labels k,X,a are 
respectively the initial momentum, initial polarization and initial color of the gauge fluctuation 
represented by a±k\ a , and the sign ± specifies whether it is a positive or negative energy wave in 
the remote past. 

2.2.2 Power counting rules 

At leading order (tree level), the energy momentum tensor is of order Q^/g 2 - In the absence of 
secular divergences, from the power counting described previously, the NLO corrections should be 
of order Q 4 S . This power counting could be obtained in eq. (12) by noting that 

a± kXa ~ 0(1), (16) 



(3 ~ 0(g), (17) 

T u ~ JJ~°(9), (18) 

since A ~ 0(g~ l ). The existence of instabilities implies that we must alter our estimate of the 
order of magnitude of the operators T„. Indeed, since T u A(t,x) is the propagator of a small 
fluctuation over the background field between a point on the initial proper time surface and the 
point (t, x), it grows at the same pace as the unstable fluctuations. Thus the counting rule for T u 
should be modified to read 

T,~0( S e^). (19) 

The combination T U T„ in eq. (12) then grows as g 2 exp(2y/JFr) which leads to a break down of the 
power counting at the proper time r max defined in eq. (7). At r max , the 1-loop correction becomes 
as large as the leading order contribution, and one may anticipate that an infinite series of higher 
loop corrections also become equally important at this time. 

2.2.3 Selection of the leading terms 

Our goal is now to collect from higher orders all the terms that are leading at the time r max . 
This comprises all the terms where the extra powers of g 2 are compensated by an equal number 
of factors of e 2 ^"\ We presume that a typical higher order correction to the energy momentum 
tensor can still be written in the form of eq. (12), but with a more general operator acting on 
T^(x) of the form 

J d 3 Ul ■ ■ ■ d 3 u n r n (ui, • • • , u n ) • T Ul • • • T Un . (20) 

Here T n is an n-point function, which may or may not be simply connected. This expression has 
not been proven in general but results from a conjecture that inclusive quantities at all loop orders 
can be expressed purely in terms of retarded propagators, thereby generalizing known results at 
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Figure 1: Representation of the 1-loop contribution involving the function Fa (it, v). The thick red 
line is the t — + surface on which the initial value problem is set up. The open circles represent 
the initial data. The filled blue circles represent the two operators T Ult) , and the U-shaped wavy 
line represented below the light-cone is the function T2(u,v). 

LO and NLO. While there are specific examples of loop contributions that have been checked to 
satisfy this conjecture, there are in particular nested loops contributions for which the conjecture 
is difficult to confirm. In the figures 1 and 2, we illustrate this formula by some examples of I-loop 
and 2-loop contributions. 

With the stated assumption implicit in eq. (20), the following power counting can be established. 
If eq. (20) is a piece of a L-loop correction to the energy-momentum tensor, the order g p and the 
number n of points of the function T n are related by 

n+p = 2L. (21) 

This formula can be checked for the examples of graphs given in the figures 1 and 2. Note that 
p = is the smallest possible value for p. Taking into account the effect of the instabilities (i.e. 
one power of exp(y / 7iT) for each of the n operators T Ui ), the order of magnitude of a contribution 
obtained from eq. (20) is 

n 

(22) 

relative to the leading order contribution. If we just count naively the powers of g, the power 
counting would indicate that this contribution gives a correction of order g 2L , a decrease by a 
factor g 2 for each extra loop. However, because of the instability, by the time r max given in 
eq. (7), the contribution is instead of order g p and does not depend anymore on the number n of 
T operators. At the time r max all the terms with p — 0, regardless of the number of loops, are of 
the same order while all the remaining terms for which p > are suppressed by additional powers 
of g. It is therefore natural to resum all the p = terms, and to neglect all those with p > as 
giving sub-leading contributions. This implies that the numbers n of T operators must be even 
and equal to 2L. Moreover, since we keep only terms of order p = in T2L, the only possibility 
that remains 8 is to construct T2L as a product of L factors T2 (an example of which is the left 

8 Note that tadpole contributions such as the term /3 ■ T u in eq. (12) are also excluded since ~ 0(g). 
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Figure 2: Representation of two examples of 2-loop contributions. The thick red line is the r = + 
surface on which the initial value problem is set up. The open circles represent the initial data. 
The filled blue circles represent operators T B ,„. Left: contribution with a that factorizes into 
two r2 ! s. Right: contribution with a F3. 

diagram of figure 2) because any non-factorized contribution to T 2 l requires more powers of g. 
2.2.4 Resummation: formal expression 

Therefore, the leading operator at L loops in eq. (20) is the L-th power of the 2-point operator 
that appears at 1-loop, 

1 L 

J d 3 u d 3 v T 2 (u,v) ■ T„T„ , (23) 
d 3 k 



1 

D. 



where 



/d 3 k 
(2?r) 3 2fc [a+kxa ■ T„][a_ fcAQ • T„] 



(24) 



The inverse 



v^erse factorial prefactor is a symmetry factor that prevents multiple counting when the 
various factors T2 are permuted. Summing all the contributions from L = (leading order) to 
L = +00, we obtain 9 



^rcsummcd( 2 ') — CX P 



i J d 3 ud 3 v r 2 (M,«)-T„T t 



(25) 



9 Until the conjecture in eq. (20) is proved, one may choose to interpret eq. (25) as a well motivated ansatz 
resulting from the exponentiation of the NLO result. 
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Eq. (25) provides an expression that resums all the leading contributions 10 at the time r = r max . 
However, this expression is very formal as it is expressed in terms of functional derivatives with 
respect to the initial conditions for the classical color fields and their time derivatives on the initial 
proper time surface E. Fortunately, as we shall see, we can rewrite this result in a form that is 
much more transparent both conceptually and for computational purposes. 



2.2.5 Resummation: path integral representation 



We first recall that the operator T u defined in eq. (13) is the generator of shifts of the initial 
conditions at all points u on the initial proper time surface r = + for the classical fields A^ 
and their time derivatives d T A^, the latter either being equal to or simply proportional to the 
corresponding electric fields, their canonical conjugate momenta. We can therefore write 



exp 



L 



d 3 u [a ■ T u ] 



(26) 



where A = (A, £) denotes collectively all the components of the initial classical field and their 
canonically conjugate momenta on the initial time surface. One should similarly understand a to 
denote small perturbations of both the initial classical field and their canonically conjugate electric 
fields. We next obtain 11 



exp 



^Jd 3 ud 3 v r 2 (u,v) ■ T U T V = J [Da] F [a] exp J d 3 u [a ■ T u ] 



with 



12 



F [a] oc exp 



^Jd 3 ud 3 v a(u) T 2 1 (u,v) a(v) 



(27) 



(28) 



In eq. (27), the functional integration [Da] is also a shorthand for integrations over all the com- 
ponents of the perturbation and of its time derivative on the initial surface. 



1() Sublcading contributions such as the p = 1 contribution g ■ gcxp(^/JW) arc no longer sub-leading by times 
r ~ ln 2 (g~ 2 ). This time is only slightly larger than the time T max at which the p = terms become important. 
Therefore, including one by one the p = 1 terms on top of an expression that resums the p = terms is bound to fail - 
these contributions must be included all at once via a resummation. Having this in mind, a more important question 
is: do these contributions ever become important after they have been appropriately rcsummed? A resummation 
of the p = 1 secular terms is outside the scope of this work, but it is plausible that they can be included in our 
framework by a modification of the distribution Fq [a] of the fluctuations at the initial time. If this is the case, then 
the p > 1 terms would simply lead to small non-gaussianities in the spectrum of fluctuations. However, the is at 
present nothing more than a conjecture, and the results of this paper have to be interpreted with this in mind. 

11 An elementary form of this identity, 



f(x) = r 



dz 



-z 2 /2a 



■ f{x + Z) 



can be proved by performing a Taylor expansion of the exponential on the left hand side and of f(x + z) on the 
right hand side of this expression. Prom this simple example, one sees that a Gaussian operator in derivatives is a 
smearing operator that convolutes the target function with a Gaussian distribution. 

12 The unwritten constant prefactor, proportional to [detl^)] -1 / 2 , is such that the distribution fbfa] has an 
integral over a normalized to unity. 
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With the identities (26) and (27) in hand, we can rewrite our formal result in eq. (25) as 

K—J*) = / [M Fo [a] T£[A + a](x) , (29) 

where the weight functional Fq, corresponding to the initial spectrum of fluctuations, is defined in 
eq. (28). This result is a central expression of our paper and was sketched previously in [18,63]. 
It was also obtained in previous works for a scalar field theory [39] and a gauge field theory [57] 
respectively using different methods. The expression is quite remarkable because it demonstrates 
that the resummation of loop (quantum) corrections that correspond to the most unstable con- 
figurations in a single heavy ion collision event can be expressed as an average over a Gaussian 
distributed ensemble of classical configurations in the Glasma. Note also that, although this for- 
mula was derived here for the energy-momentum tensor, the power counting that led us to the 
exponentiation of the 1-loop result did not depend on the choice of a specific observable. Thus 
we expect that the same resummation would also be applicable to other inclusive quantities; the 
spectrum of fluctuations superposed on the Glasma fields is universal. 

The essential ingredient in eq. (29) is the small fluctuations correlator ^(m, v), defined in 
eq. (24), which should be computed with the two endpoints on the initial time surface S, with the 
Glasma field in the background. A first attempt to compute this object is given in [57]. However, 
the expression obtained there was incomplete because it corresponds to an approximate expression 
for the free propagator in A T = gauge, with the only dependence on the background field coming 
from Gauss's law. We will compute here (in A T = gauge) the spectrum of fluctuations both in 
the free case and in the presence of a classical background field. We will show later that the latter 
has a non-trivial dependence on the classical fields in the Glasma at E. 



2.3 Gauge invariance of the spectrum of fluctuations 

The left hand side of the expression in eq. (29) should be gauge invariant because the energy- 
momentum tensor is a physical quantity. It should therefore be invariant under a gauge transfor- 
mation of the classical Glasma field, 

A — > nU**n + ^ft^ft . (30) 

This invariance is also true for its leading order counterpart on the right hand side of the expression, 
when we apply the same gauge transformation to the total field A + a, 

A + a^rt (4" + a") ft + ^nt0"ft . (31) 

Thus for the expression in eq. (29) to be manifestly gauge invariant, it is sufficient if the initial 
spectrum of fluctuations F is invariant under the transformation 

a — > ft f aft . (32) 

If we decompose a on the basis of the generators of SU(3), a = a a t a , the identity 

n ab t b = ft t a ft t , (33) 
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(where f2 a f, here is an adjoint SU(3) matrix) gives the equivalent transformation for the components 
a a to be 

a a -> Cl ab a b = & a . (34) 
For the argument of the exponential in F to be invariant under this transformation, 

a a (u) T~^ ab (u, v) a b (v) — > a a (u) T^ ah {u, v) a b (v) , (35) 

the inverse small fluctuations correlator in the Glasma (which is an 8x8 adjoint matrix in SU(3)) 
must satisfy 

T^(u,v) =n(u)T2 1 (u,v)n jt (v) . (36) 
It is clear from the structure of the expression in eq. (24) that this property will be satisfied. 

2.4 Renormalization of ultraviolet divergences in the Glasma 

It is also important to address the potential ultraviolet divergences in cq. (29). The leading order 
energy-momentum tensor in the Glasma is ultraviolet finite. However, because one is resumming 
quantum fluctuations in eq. (29), the energy-momentum tensor should receive a contribution from 
the (infinite) zero point energy. One can regularize ultraviolet divergences by introducing a cutoff A 
corresponding to the largest momentum mode 13 of the fluctuation a. Since the energy-momentum 
tensor has canonical dimension four, we expect that its dependence on this cutoff can be organized 
as 

^ c r_ d W- C iA 4 + c 2 A 2 + c3, (37) 

where 01,2,3 are finite quantities. It is easy to renormalize the energy-momentum tensor by subtrac- 
tion if one can prove that the divergences are truly a property of the vacuum and do not depend on 
the background classical field A in the Glasma. The coefficient c\ is dimensionless - it is therefore 
a pure number, that cannot depend on the background field. The case of c 2 is trickier. Indeed, its 
canonical dimension 2 allows a priori a dependence on the background field. However, we know 
that the left hand side in eq. (37) is invariant under gauge transformations of the background field; 
we therefore must conclude that the coefficient ci must be a gauge invariant, local (because the 
left hand side is a local quantity), dimension 2 quantity. There is no such quantity in Yang-Mills 
theory, which suggests that C2 = 0. Thus, on the basis of gauge symmetry and locality, we expect 
that the only ultraviolet divergence in our expression for the resummed energy-momentum tensor 
is a quartic divergence, with a coefficient that does not depend on the background field. It can be 
computed in principle once and for all in the absence of the background field and subtracted from 
^resummed to §i ye an ultraviolet finite result for this quantity. 

3 Orthonormal basis of small fluctuations 
3.1 Introduction 

We noted in eq. (29) of the previous section that the small fluctuations 2-point correlator r 2 (it, v) 
is the key ingredient in resumming the contributions of leading instabilities at all loop orders to 
the energy-momentum tensor. 

13 In practical implementations of this resummation, space is discretized on a lattice, and thus the UV cutoff is 
the inverse of the lattice spacing. 
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In this section and the following one, we will work in the Fock-Schwingcr gauge A T = 0. 
Albeit at first sight a natural gauge for describing hadron-hadron collisions (because it is an 
interpolation between two light-cone gauges in the forward light-cone) , the Fock-Schwinger gauge 
is not frequently used in the literature. This is because even expressions for the free correlator are 
complicated in this gauge. Our motivation here is specific to the nature of the CGC description of 
heavy ion collisions. The initial conditions for the evolution of classical gauge fields in the forward 
light cone, in this gauge, are simply expressed [47] in terms of the classical fields in the nuclei before 
the collision. This is an important criterion because our resummed result in eq. (29) for the energy- 
momentum tensor is expressed in terms of solutions of classical Yang-Mills equations with Gaussian 
distributed initial conditions. Further, numerical computations are unavoidable because one is in 
a strong field regime where perturbative computations are invalid; thus analytically cumbersome 
expressions are not a deterrent if efficient numerical algorithms are feasible. 

Turning to the computation of the small fluctuations correlator, 

/d 3 k 
( 27r )3 2 fc a +k xa{u)a- kX a(v) , (38) 

as noted previously in the discussion after eq. (13), the fluctuation fields a^ feAa are plane wave 
fields at x = — oo that have been evolved in time by interacting with the classical background 
field A. Before going further, let us state two properties of this correlator that are true when the 
two points u and v lie on the same Cauchy surface 14 , 

i. T 2 is symmetric: r 2 (it,t>) = r 2 (f,M), 

ii. r 2 is real valued. 

Note that for the correlator as defined, the two properties are equivalent since ci-kXa = ( a +kXa)* ■ 
These properties are crucial since T 2 is the variance of the fluctuations in our resummation formula, 
eq. (29). 

From the definition given in eq. (38), one might presume that the calculation of T 2 requires one 
to follow the entire evolution from x° = — oo to r = + . We will demonstrate in section 3.3 that 
this is not the case; to compute T 2 , it is sufficient to construct an orthonormal basis of fluctuation 
fields about classical fields at small proper times r = + . 



3.2 Curvilinear coordinates 

The fluctuations a± kXa start as plane waves in the remote past and evolve about the Glasma 
classical field A according to the equation of motion (14). This expression of the equation of 
motion, as well as the Yang-Mills equation (11), implicitly assume a Cartesian system of co- 
ordinates. However, the most natural co-ordinate system in the treatment of the post-collision 
evolution is the (r,r],x±) system, where 

r = Vt 2 - z 2 , V=l^(^^j , (39) 

14 For these properties to hold true, it is important that the two points are not separated by a time-like interval, 
which ipso facto is guaranteed if they belong to the same locally space-like surface. 
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and x± collectively denotes the two co-ordinates perpendicular to the collision axis. In these 
co-ordinates, the metric tensor 

^ = diag(l,-l,-l,-r 2 ) (40) 

has a proper time dependent determinant, y/— g = t. This implies small changes to the equations 
of motion. The classical Yang-Mills equations (11) become 15 

~^=f a W~gg a0 g vli T^) = S v + Pl + 5 V ~ p 2 , (41) 

and the small fluctuation eqs. (14) become 

~^=D a (J—gg af3 g^ (V^ - V^)) - igg a0 g^^ a a =Q. (42) 

In Fock-Schwinger gauge A T = and (r, r], x±) co-ordinates, these equations can be written out 
more explicitly to read, 

V v d T a v + T 2 Vid T ai - ^d T T> v a, q - ir 2 a T 2? i a i = 

(d T T~ 1 d T - t^ 1 Vo) a v + T~ l V iv ai = 
(<9 t t<9 t - T~ x V m - tVu) a x + t~ V vx <h\ + rV lx ai = 
(d T Td T - T^Vnr, - tVu) a y + r'^V^a^ + TV iy ai = 0, (43) 

where we have introduced the shorthand V, which is defined as 

Kj< = V a ;vfa b K + gf acb ^ ja b K . (44) 

Here the indices I,J 7 K denote x,y and 77, while we use Latin indices i,j to designate the two 
transverse coordinates x,y. For example, d 2 — d 2 + d 2 = d\. Note that the first equation in 
eq. (43) is Gauss' law, which can also be written as 

V, q d T a, q + T 2 V,d T a, + ig (d T A v ) a. q + t 2 {d T A t ) a t = . (45) 



lj Field equations can be generalized to an arbitrary system of coordinates by trading ordinary derivatives 9 M for 
covariant derivatives V M (here, covariant refers to coordinate transformations, not SU(3) gauge transformations) 
that involve Christoffel symbols. For the Yang-Mills equations, one should thus write 

(V „ - igA^T^ = J v , = VM" - VA" -ig[A»,A v ] . 

However, it turns out that for an antisymmetric tensor such as one has also (see [64], chapter 5) 

1 



v — 9 

Moreover, one can show that V^A,, — V„A M = c9 M A„ — <9„A M , so that the flat-space expression of J-" M „ can still 
be used (provided the two indices are downstairs). In other words, the usual formulas for the field strength should 
be used only for lower indices, and the metric tensor should be used to raise indices if necessary. For instance, in 
A T = A T = gauge, one has T Tr] = d T A v and T rT > = g Tr g^ T TTI = -r- 2 d T A v , but ST A 7 ! = d T (-T- 2 A v ) is not 
equal to T TTl ' . 
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Even though we exclusively work in r-r) coordinates, we shall at times make use of light-cone 
coordinates defined as 

(46) 



x ± = 



V2 



with the metric tensor now having the form 

.9+- = .9- + = 1 ! 9xx = 9yy = -1 • 

The transformation from light-cone to r, r\ coordinates is given by 



re 



±7, 



X = 



V2 ' 



while the Fock-Schwinger gauge condition, expressed in light-cone coordinates, is 

A T = A T = - (x+A- + x~A+) = . 



(47) 



(48) 



(49) 



3.3 Inner product for small field fluctuations 

Since the equations of motion (42) of the small fluctuations are linear, the set of its solutions is 
a vector space, and it is sufficient to know a basis of this space in order to be able to construct 
any solution. For a real background field such as the classical field >4 M , the evolution in time of 
the small fluctuations is unitary 16 . Therefore, there should be an inner product between pairs of 
solutions of eq. (42) that remains invariant during the evolution of these solutions. To construct 
this inner product, rewrite eq. (42) as 

O^a^ = , (50) 

(51) 

d 4 xa* u (x) O v »-O v »* b^x) , (52) 



with 



Consider now two solutions a M and 6 M of eq. (50), and start from the identity 







/ 



where Q is some domain of space-time. This identity is a trivial consequence of the equations of 
motion for a* and b. A remarkable property of the integrand in the right hand side is that it is a 
total derivative 17 , 



a*(x) 



Ql>H _ Qh 



9\9 vll 9^ 



Therefore, one can use Stokes theorem, 

d A x d n F a 



g ^g^- \^ a ^)[al{x)%b^x))] . (53) 



f d 4 x d a F a = [ , 
Jn Jan 



d 3 S u n a F a 



(54) 



16 In event of confusion from the apparent structure of the last term of eq. (42), note that the components of the 
adjoint generators are purely imaginary, and therefore the function that multiplies a M in this term is real. 

17 For this property to be true, it is crucial that the last term in eq. (51) is real and one should properly take the 
complex conjugate of the covariant derivatives when they act on the left. This property is in fact closely related to 
the operator 0" M being Hermitean; the evolution of the fluctuations is unitary. 
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where d 3 S u is the measure on the boundary dCt, and n a is a normal vector to the boundary, 
oriented outwards. Let us assume that the boundary <9f2 is made of two locally space-like surfaces 
Si and E2, and a third boundary located at infinity in the spatial directions on which all the fields 
are vanishing. Then eq. (52) is equivalent to 

d 3 s u V=g(g^g aP - \<? p <r - \g va 9^)n a (a* v {u) v fi &») 

= / d^S u ^ g [g^g a ' 3 -\g^g^ a -\g m g^)n a (al{u)Vpb^u)) . (55) 

We have thus proved, most generally, that an inner product defined as 

(a\b) = ij <?S U V—g(g^g a ' 3 - ^<T - \g va g^)n a (a* v {u) V b^uj) , (56) 

is independent of the Cauchy surface S used to define it, provided a p and b^ obey the equation of 
motion of small fluctuations. Note that we have added a factor i to its definition to ensure that it 
is Hermitean, 

(a 1 6)* = (b\a) , 

(a*\b*) = -(b\a) = -(a\b)* . (57) 

In the special case where S is a surface of constant r and we work in the Fock-Schwinger gauge 
A T = 0, we have n-V = d T , and n ■ a = 0, n ■ b = 0. Therefore the inner product simplifies into 

(a\b) = i J d 3 s u V^gg^(<(u)d T ■ ( 58 ) 

r=const 

Now let us evaluate the inner product for pairs of field fluctuations taken from the set of the 
a±k\a- Since the inner product does not depend on the chosen time surface and since we know 
these fields at x° — > —00 (because they are defined via their initial condition in the remote past), 
we can evaluate the inner product by using plane wave initial conditions for these fluctuation fields. 
This gives 

(a+k\a\a-i p b) = 

(a+k\a\a+i P b) = 5 Xp 5 ab (2TT) 3 2kS(k - 1) 

(a-k\a\a-i pb ) = -5 Xp 5 ab (2ir) 3 2k6(k - 1) . (59) 

Thus this particular basis of the space of solutions of eq. (14) is orthonormal with respect to the 
invariant inner product defined in eq. (56). Note also that the a+k\a's represent only one half of 
the basis of the vector space of solutions of eq. (42) -namely the solutions that have a positive 
frequency in the remote past. The other half is simply obtained by complex conjugation. It easy to 
check that any unitary transformation of the positive energy solutions (and a concomitant change 
to the negative energy ones, that are their complex conjugates) transforms an orthonormal basis 
into another orthonormal basis, and leaves the formula eq. (38) unchanged. This remark is useful 
because it leaves us the freedom to label the elements of the basis by other quantities than the 
Cartesian 3-momentum. This will be true in our specific case where we are interested in a basis in 
a curvilinear co-ordinate system. 
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3.3.1 Normalization of the fields and choice of basis 



It is important to note that the prefactor in front of the S functions in eq. (59) exactly cancels the 
factors that are included in the integration measure in eq. (38), namely one has 

d 3 k 

(27r)32fc ( a + fcA °l Q +'P b ) = 1 ■ ( 6 °) 

This remark in fact defines uniquely how the inner product of the basis elements should be normal- 
ized given a generic choice for the integration measure in eq. (38); in particular, this rule will be 
come in handy later when we use other labels than the usual 3-momentum to index the elements 
of the basis. Moreover, this makes clear that eq. (38) is just one particular representation of the 
correlator r 2 ; there exists such a representation for any orthonormal basis of the space of solutions 
of eq. (42), as we shall explain now. Thanks to the above inner product, one can spell out a general 
procedure 18 for constructing the correlator T 2 : 

i. Find a complete set of independent positive energy solutions a K of eq. (42), where K denotes 
collectively (usually a mix of continuous and discrete labels) all the labels necessary to index 
these solutions. 

ii. This set of solutions should obey the orthogonality condition, 

{a K \a K ,)=N K 8 KKl (61) 
with N K real and positive definite 19 , 
iii. The correlator r 2 is then given by 

r 2 (u,v) = J d f i K a K (u)a* K (v) , (62) 

where the measure dfi K (a mix of integrals and discrete sums) is such that 

d» K N K S KK , = 1 . (63) 



It is clear from eqs. (61) and (63) that the r 2 given by eq. (62) is independent of how we normalize 
the solutions (i.e. on the constants N K ), provided we choose the integration measure accordingly. 
Moreover, we only need to know the form of the solutions at the time of interest, and we can avoid 
the complication of evolving the plane waves from the past through the forward light-cone. This 
is particularly helpful when one uses r, r\ coordinates, because this system of coordinates has a 
singularity at r = 0. 

A further simplification is possible because in practice we won't need to use directly eq. (62) 
for T 2 . Indeed, an ensemble of real- valued field fluctuations a M that have a 2-point equal-time 
correlation given by T 2 can be generated by the following formula, 

a»(x) = [ d^ K \c K < (x) + c* <*(*)! , (64) 



18 In this light, eq. (38) which represents T2 in terms of the a±kx a 's, exploits one possible method of constructing 
such an orthonormal basis. In this case, one starts at x° = — oo with the plane waves, that are known to form an 
orthonormal basis, and evolves them forward to the time of interest. The time invariance of the inner product then 
guarantees us that we get an orthonormal basis on the forward light-cone. 

19 This means that the solutions a K will in general be complex solutions. 



18 



where the coefficients c K are random Gaussian-distributed complex numbers whose variance is 
given by 

(^ K ) - ^s KK , 

= («,)=0. (65) 

This method of generating the field fluctuations offers the advantage that it does not require that 
one diagonalizes the correlation function T 2 ■ 



3.4 Free fluctuations in A T = gauge 

Let us start by calculating the fluctuation correlator in A T = gauge and on the surface r = + 
in the free case, namely, in the absence of the classical background field. Given the complications 
introduced by the choice of gauge and the system of curvilinear coordinates, this is a useful exercise 
to pursue before attacking the more general case of the Glasma background field. In this situation, 
eqs. (43) simplify into 

d v d T a 71 + r 2 did T ai = 
(d T T~ 1 d T — r~ 1 d 2 1 ) ciri + T~ 1 did, 1 a i = 
(d T Td T - T~ 1 d 2 - Td\) a x + r~ 1 d n d x a n + rdid x ai = 
(d T Td T - T~ 1 d 2 - rd^) a y + T~ 1 d, 1 d y a v + rdidydi = 0. (66) 



3.4.1 Residual gauge freedom 

A general solution a M to eq. (66) has a priori four components. However, a massless vector field has 
only two physical degrees of freedom. One of the seemingly independent components of the vector 
field is removed by the gauge condition a T = a T = (this is already implemented in eqs. (66)). 
However, even after imposing this condition, there is a residual gauge symmetry in the equations 
of motion, namely, these are invariant under r independent gauge transformations 

a" -> +d ll A(r 1 ,x ± ) , (67) 

where A is an arbitrary r independent function. As a consequence of this residual gauge freedom, 
the three remaining components of a M are not all physical degrees of freedom. 

In order to find the two physical solutions, we begin by finding the unphysical solution. This 
solution must be a pure gauge, which here means it is a r independent total derivative. As we 
will see, the unphysical solution is not a dynamical variable but is completely constrained by the 
initial and boundary conditions. After finding the most general r independent solution to the 
equations of motion, the two physical solutions can be determined relatively easily. Their form 
will be narrowed down by requiring that the three solutions are mutually orthonormal, and then 
the residual gauge freedom will be fixed by imposing the equations of motion. 

In the (r, r), x±) system of coordinates, a convenient set of labels for the solutions of eq. (66) 
is v, fcj_, A, a, where v is the Fourier conjugate of the space-time rapidity r\ (as used previously, A 
denotes the polarization and a the color). We choose A = 1,2 to denote the physical solutions, 
and A = 3 to be the unphysical one. Since in this section we are in the vacuum, all the colors a 
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obey exactly the same equations and we can drop this index to keep the notations lighter. With 
this set of variables, eqs. (61), (62) and (63) take the form 

T 2 (x,x') = E /^au(r = 0+xK A (r = 0+x'), 

(akA|ak'v) = (2ir) 3 6„ x , 8(v-v')8 ( k ± -k' 1 _) , (68) 

^ V " ^^^^^^ 

a(k-k') 

where we use the shorthands k = (v, kj_) and x = (77, x_\_). 
3.4.2 Vacuum solutions 

In a linear system of coordinates, the a^ x (x) introduced above would have the simple following 
parametrization, 

<aW=^ W 4, (69) 

where e^ x is a constant polarization vector. Note that the minus sign in the exponential that gives 
the time dependence is necessary in order to ensure that a^ A is a positive energy solution. But 
because we work in a curvilinear co-ordinate system, the time dependence of the solutions cannot 
be a simple exponential. Let us generalize the previous expression by writing 

< A (T,x)= e * k - X < A (T), (70) 

where we have combined the polarization vector and the time dependence in a unique quantity that 
we denote a^ A (r). (Since the equations of motion do not have coefficients that depend explicitly 
on rj or x± 7 it is clear that we can still look for solutions whose 77 and x± dependence is of the 
form exp(«k • x).) Here again, we will have to make sure that the a^ A (r, x) constructed in this way 
contains only positive energy contributions. 

Let us consider first the unphysical solution. This is a pure gauge solution independent of r. 
The most general r independent solution is of the form 



*k3 



k y I <!;,(;/. fc ± ) . (71) 



V 



where a^(y, fej_) is an arbitrary function. The inner product of the unphysical fluctuation a£ 3 with 
one of the physical solutions (a£, A with A = 1, 2) is 

/+00 ,■ 
dt] / d 2 x ± e «( k '- k )- x d T (k x a£, x + k y al, x + vT- 2 al lx ) 
-OO J 

= zr(27r) 3 (5(k - k')a* 3 (u, k ± ) 8 T (k x a£ x + k y a y kx + ^- 2 a£ A ) . (72) 
We can satisfy this orthogonality condition by choosing a second solution of the form 

k 

<j(r) = I -L j n ir u.k_ ) . (73; 
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The functional form of a\ should be fixed by the equations of motion. Substituting the above 
expression into the equations of motion (66) yields a differential equation for a±, whose general 
solution can be expressed in terms of the Hankel functions and 

ai(r, u, k±) = a k (fcjr) + 6 k (k ± r) . (74) 
Recall here the integral representation of the Hankel functions, 

H[l\x) = e +7rl//2 / e «cosht+^t^ 



/-\-oc 
e -ixco S ht-iut dt ^ (75) 
-oo 
(2) 

The convention set by eq. (69) implies that only H- v '(k±T) has the appropriate frequency sign to 
be one the akA's. We can therefore set a k = and keep only the second term in eq. (74). The 
value of 6k can then be determined by the orthogonality condition. For this, we need the identity 



from which we obtain 



H^( X )d x Hl?( X ) = -^, (76) 
(a k i|a k a) - (2^) 3 <5(k - h')\b k \ 2 kl^-^ . (77) 



In order to get the same normalization as in eq. (68), we obtain (up to an irrelevant phase) 



W = ^— , (78) 



and the first physical solution reads 



<i(x) = [ X ] e*-ffW(fc ± r) . (79) 







The second physical solution can be found by requiring that it be orthogonal to the two solutions 
we have so far. This restricts it to be of the form 

(vk x a 2 ±(r, v,k±) \ 
vk v a 2 ^{T,v 1 k 1 _) , (80) 
-a 2v (T,v,k±) J 

provided that k 2 _T 2 d T a2± = d T a 2r) . This last constraint was derived by substituting the general 
form of the second solution eq. (80) into the orthogonality condition cq. (72) with the unphysical 
solution. It turns out that this is the same constraint needed for this solution to fulfill Gauss's law. 

Dynamical equations for the functions a 2 ± and a 2r) can be found by substituting eq. (80) into 
cqs. (66). One obtains, 

d z T a 2ri - ^d 2 a 2n + + k i^j d T a 2n = , 

d*a 2± + ^d 2 a 2± + + fel) d T a 21 _ = . (81) 
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The positive energy solutions to the above third order differential equations can be written as 



a 2r) = const x J dr'r' (k±T r ) 

a 21 _ = const x J* ^-H^(k^r') . (82) 

Note that the differential equations above imply that the functional form of a 2 ± and a 2v are 
only determined up to an arbitrary r independent function. This corresponds to a residual gauge 
freedom in which we can always add to the second physical solution cq. (80) a pure gauge solution 
having the form of eq. (71). This residual gauge freedom can be removed by imposing an additional 
gauge fixing condition. For example, in simulations of classical Yang-Mills one typically imposes 
transverse Coulomb gauge. 

The integrals over the Hankel functions can be written in terms of hypergeometric functions 
but their form is not very enlightening. To streamline our notation, following Makhlin [58], let us 
define 

= f dx x b H<?\k ± x) . (83) 
The properly normalized form for the second physical degree of freedom is 

a^{r,r,,x ± ) = I uk v R% v (fc ± r) | e- x . IH4 ; 



For k±r <C 1, we can make use of the series expansion of the Hankel functions and rewrite the 
solution as 



k x 
k 

2k± \ -{k 1 _Tf/{y + 2i) 



<2««4f^( k y ) e*-* Hg\k ±T ) . (85) 



From these explicit solutions, we will construct in appendix B the correlation function T 2 for free 
fields on the initial surface r = + . However, for the purposes of generating a Gaussian ensemble of 
fluctuations with the proper variance, the above results along with eqs. (64) and (65) are sufficient. 



4 Small fluctuations in the Glasma 

After our extended discussion of the free fluctuations in A T = gauge, we now turn to the derivation 
of the small fluctuations spectrum in the Glasma. We shall first write down the small fluctuation 
equations of motion in the presence of the background classical fields. The main difficulty here is 
that the Glasma background fields are not known analytically at arbitrary proper times. They are 
however known in closed form at t = + in terms of the classical CGC fields before the collision. 
We will perform a small time expansion, valid at very short proper times r <C Qj 1 , of both the 
classical fields and the small fluctuation fields and show that the fluctuations only depend on the 
classical gauge fields immediately after the collision. As our final result, we will obtain explicit 
expressions for an orthonormal basis of small fluctuations that generalize eqs. (79) and (84) to the 
case of a non-zero background field at small proper times. 
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4.1 Structure of the Glasma background field 

In the Fock-Schwinger gauge, the classical gauge field configurations can be expressed as [47,51] 

A i = e(-x + )e(x-) a \ + e(x + )e(-x-)ai + e(x + )e(x-)A l 

A' n = 9(x+)9(x-)A n (86) 

The fields a\ 2 are the color fields of the two nuclei before the collision, that take the form of 
transverse pure gauge fields, while A^ denotes the classical fields after the collision. 

Since we are interested in the spectrum of fluctuations at r = + we need only the behavior of 
the background fields shortly after the collision. The classical Glasma fields in the forward light 
cone can be expanded, in all generality, at early times as 

oo 

A T = J2 A (n)ir n (87) 

n=0 

The initial conditions for these background fields at r = + are obtained by matching the Yang- 
Mills equations just below and just above the forward light-cone (to ensure a regular behavior of 
the field equations). One obtains [47] for the fields and their time derivatives, 

A{t = 0+) = a\+ai 
A,(r = 0+) = 

£i(T = 0+) = Td T A\r=0 = 

£"(t = 0+) = \A v \ T=0 = tg[a\,a 2 ] . (88) 

As alluded to previously, explicit analytical solutions are known for the fields a\ 2 - The Taylor 
expansions of eqs. (87) begin with 

Ai=Ai(0 + ) + O(r) , ^ = ir'(0+)r 2 +O(r 3 ). (89) 

At this point it is useful to introduce some extra notation. Based on the Taylor expansion in 
eq. (87) of the classical background field, it will be convenient to introduce an analogous expansion 
of the covariant derivatives 20 and projectors as defined in eq. (44), 

n 

?v = E r,l7 w- ( 9 °) 

n 

Naturally, the coefficients T^( n )n an d ~P( n )^ v can be written in terms of the Taylor coefficients of the 
classical background field, A[ n )i- For example, if we perform the Taylor expansion of the covariant 
derivative, I?(„) AI , the leading coefficients can be expressed in terms of the Taylor coefficients of the 
classical background field, 

V#)i = 8 ab d i -igA# )i , V$ )ri = 5 a % , V^^-igA^. (91) 



J Of course, the D T derivative does not need to be expanded, since in the Fock-Schwingcr gauge one has D T = d T 
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In writing down the above expressions we have used the fact that the transverse components of the 
background field have non-vanishing zero order Taylor coefficients in contrast to the longitudinal 
component of the background field whose leading behavior starts with r 2 . Similarly we can express 
the Taylor coefficients of Vi n )nv in terms of the Taylor coefficients of the classical background field. 
The terms which will needed later on in our discussion include, 

1 (0)r?i - 1 (O)mj — U {i))i°^ ' 

V[m = Sa " d ^ - - igAfajdi - ig (d l A m ) ab g'A^A^ , 

Ppo™ = -KgAfandr, ■ (92) 

In deriving the above expressions we have used the fact that the background field is boost invariant, 
dr/Ai = 0. It is this property of the background field that leads to the simple form of V(o)rii and 
yields the factor of two in the expression 21 for V(2) m - 



4.2 Early time behavior 

In the case of the r, r], x± system of coordinates and in the Fock-Schwinger gauge, the equations 
of motion for the fluctuations propagating over such a background field are written explicitly in 
eqs. (43). Solving for the full time dependence of the fluctuations is both intractable (since the 
time dependence of the background field is only known numerically) and unnecessary (since we 
only need the spectrum of fluctuations at early times). 

A crucial property of the background Glasma fields is that they are invariant under boosts 
in the longitudinal direction. This implies that the fields Ai and A- n in the forward light-cone, 
after the collision, are independent of the space-time rapidity r\. Thus the variable v, the Fourier 
conjugate of 77, is a conserved quantum number for fluctuations propagating over the Glasma fields, 
that can be used to label the elements of the basis. We can therefore write the elements of the 
basis as 

a^ x (T, V ,x ± ) = e^^ lx (T,x ± ), (93) 

where A = 1, 2 is a polarization index and I collectively represents the remaining quantum numbers 
necessary to label the basis. The main difference with the free case (see eq. (70)) is that we cannot 
assume that the x± dependence of the fluctuations has the form of a plane wave and instead 
represent it more generally as the function /3^(t, x±). This is because the background field has a 
non-trivial dependence on xj_. A further consequence is that, in contrast to the vacuum case, the 
remaining quantum numbers encoded in I will not simply be transverse momenta. As in the free 
case, we shall keep only solutions that have positive frequencies in this basis. 

We would now like to motivate one of the main results of this work; a modified form of the 
linearized equation of motion eq. (43), which captures the early time behavior of a quantum 
fluctuation propagating on top of the classical background field. The simplest way to arrive at our 
result is to simply replace all of the projectors appearing in eq. (43) with the corresponding zeroth 
order Taylor coefficient from the proper time expansion given in eq. (90). Following this procedure 
gives essentially the right result except for one subtlety. We will argue that we also need to include 

21 It may be instructive to work this term out explicitly. Since the field strength tensor is anti-symmetric we have 

V$ = V^Vf. We can therefore write V$ )m = ^S) V V (2) V + V (2)^(0) v From °^ ( 91 ) wc know that V (0)r, = 

6 ab d n and therefore commutes with the boost invariant background field. The final result is then P" 6 , = 2X>?A c9„. 
' (2)1717 (2 V 'I 
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one term, t^V^)^, that appears at higher order in the expansion of the projector. The resulting 
small-time linearized equations of motion are, 

(d T T^ 1 d T - r _1 P( )ii) a v + r _1 'P( 0) ^a i = , 

(<9 T T<9 r - T~ V{Q) m - T V(2) m - T P(0)ii) a x + T~ 1 V(0) rix a v + T V^ix^i = , 

(d T Td T - t~ 1 7 7 (o) w - T~V(2) m - tV( )u) a v + T~ x V(Q) m a v + rV^iyO-i = 0. (94) 

The necessity of including the term {V(2) m ) can be seen by looking at the structure of the operator 
acting on a x in the second equation (or equivalently acting on a y in the third equation) above; 

(d T Td T - T~ ll P(Q) m - rV(2) m - tV( )u) a x , (95) 

By examining the power counting in r of each term in this operator, it would be inconsistent to 
ignore the t 2 component of P vv which is of the same order as V(o)u- One could argue that both 
V(o)a and Vin\ m are suppressed by r 2 at early times relative to ViQ) m and therefore can both be 
ignored. If we drop these terms, when we turn off the background field, we would not recover the 
vacuum wavefunctions and this is clearly unsatisfactory. We therefore conclude that we need to 
include all projectors that are of the same order in r as the constants appearing in the vacuum 
case. By the argument presented we should also include the term V^^-q m order to have the 
correct power counting. 

A more formal way of arriving at eq. (94) is by considering the series expansion of the small 
fluctuations, aj. Using the method of Frobenius one finds that the leading r behavior of the 
transverse components behaves as a* ~ t iu while that of the longitudinal component goes as 
a n <~ t 2+w . These coefficients are exactly those needed to reproduce the small r expansion of the 
physical solutions found in the vacuum case. Furthermore, if we continue to use the Frobenius 
method we find that only a small subset of the Taylor Coefficients in eq. (90) are needed to 
determine the lowest order coefficient in the Frobenius expansion of aj. These are precisely the 
Taylor coefficients that have been included in eq. (94). 

Finally, we need to discuss why we have not included Gauss's law in eq. (94). The reasoning 
is that Gauss's law is not a dynamical equation but a constraint, and therefore not amenable to a 
series expansion. This can be seen by noting that Gauss's law, 

Q = t- x V^ + tD^ 1 = (96) 

is a constant of motion {d T Q — 0) and therefore if Gauss's law is obeyed at t = + it will remain 
satisfied for all times. 



4.3 First physical solution 

For the first physical solution, we take = as was done in the vacuum case. With this choice, 
the first equation in (94) coincides with Gauss's law. The last two equations control the evolution 
of the transverse components a l vll . Since we expect the time dependence in eq. (93) to enter in the 
same way as in the vacuum case, we postulate that the solution at early times will be of the form 

aZ n (T,r ] ,x ± )=V— bl n { Xl _) e™> Hg> (Q vll r) . (97) 
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Next we substitute eq. (97) into the early time linearized equation of motion (94). Requiring that 
cq. (97) is a solution to the equations of motion leads to, 



V {o)y V (o)y + 



(") 

(2)w 
.(") 

(2)W 



ull ■ 



(98) 



We have explicitly included a superscript (v) on V" 2 ) m m order to remind ourselves that the 
derivative with respect to rj should be replaced by iv when acting on the exponential in eq. (97). 
We should also stress that the equations above only depend on the background fields at r = + , 
which are known analytically. 

Solving eqs. (98) amounts to finding the eigenvalues Q 2 n and eigenfunctions (the doublets 
(bf, ll (x±),b 1 ^ ll (x±))) of an Hcrmitcan operator. Since this operator is Hermitean, its spectrum is 
made of real eigenvalues, and its eigenfunctions can be chosen to form an orthonormal basis, 



/ 



d 2 x ± bi* n (x ± )bi vi {x ± ) = 6 U > 



(99) 



Note that the choice of normalization in eqs. (97) and (99) is such that the inner product defined 
in eq. (68) is satisfied 

{a vll \a v , vl ) =2n5(u-i/)5 u > ■ (100) 

The operator that we need to diagonalize in eqs. (98) has a spectrum that is twice larger than 
the size expected for the space of solutions with polarization A = 1. Half of this spectrum is 
incompatible with Gauss' law and must be discarded 22 . 



4.4 Second physical solution 



As was the case for the first physical solution, the second physical solution will maintain the same 
t dependence at t = + but with a modified dispersion relation. We therefore write the most 
general form of the second physical solution as 



"vl2 



TTv/2 



2Q 



vl2 



( b- l2 (x ± )R% v (Q^r) \ 

bl l2 (x ± )R% v (Q vl2 r) 
V bl l2 { Xl _)R%„(Q ul2 T) ) 



(101) 



Following the same procedure as for the first solution, we substitute eq. (101) into the linearized 
equations of motion given by eq. (94). Requiring that eq. (101) is a solution to the equation of 
motion at lowest order in r leads to the following equations for b l ul2 



12 



12 



wV ( Q )x bl l2 
ivV (0)yKl2 



(102) 



22 In the free case, this operator is Oij = —d\&ij + didj. It has two types of eigenfunctions: (i) b' = diXi with 
eigenvalue Q 2 = 0, and (ii) b 1 = tijdjX, with eigenvalue Q 2 = k'j_ (where x( x ±) = cxp(ik± ■ x±)). Only the second 
eigenfunction is compatible with Gauss' law dib 1 = 0. Interestingly, one of the perturbative solutions of the classical 
equations of motion at large transverse momenta in the forward light-cone has an identical structure [47]. 
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The first two of these equations simply give in terms of &^ i2 - The third equation determines 
b v ul2 as an eigenfunction of the operator —V(p)u, with eigenvalue Q^ 2 - Because this operator is 
Hermitean, these eigenvalues are real, and the eigenfunctions are mutually orthogonal, 



(103) 



We can therefore write the final form of the second physical solution in terms of the single 
eigenfunction b n ul2 



,W2 



2Q 



vl2 



V 



(2) 

ivR -l.iu {Qvl2T) V {0)y 



Kl2(x±V V - 



(104) 



where ^ i2 (j:i) is a solution to the eigenvalue equation (102). Finally, let us rewrite the solution 
using the small time approximation of R±\ iv , as done in the vacuum case (see eq. (85)) 



< 2 (x) 



v/2 



2Q„ 



12 



v {0)y 

-(Qvl2T) 2 /{v + 2l) 



(105) 



5 Outline of an algorithm for numerical computations 

The results of the previous section provide all the ingredients we need in order to evaluate an 
inclusive quantity such as the energy-momentum tensor, resumming in the calculation both the 
large logs of l/xi,2 and the secular terms that plague fixed order calculations. The algorithm for 
performing such a calculation can be broken down into several independent steps: 

i. Solve the JIMWLK equation. 

i.a Generate an ensemble of color source densities p a (x±) (or, equivalently, of Wilson lines 
Q a b(x±)) that represent the distribution W[p] at large x, close to the fragmentation 
region of a nucleus. 

i.b For each of these configurations, evolve it to smaller x by using the Langevin formulation 
of the JIMWLK equation. This amounts to performing a random walk on the space of 
mappings from R 2 to the group SU(3) [65,66]. 

ii. Pick two elements (one for each nucleus) in the above ensembles for each projectile, evolved 
at the values of x relevant to the observable of interest. Compute the gauge fields and their 
first time derivatives on the initial surface t = + immediately after the collision [52,56]. 

iii. Generate fluctuations on top of the classical Glasma fields. 

iii.a Solve the eigenvalue equations in eqs. (98) and (102). In the former case, only the 
solutions that fulfill Gauss's law should be kept. It should be noted that in a lattice 
discretization, this amounts to diagonalizing large but sparse matrices. 
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iii.b Evaluate the Hankel functions H\j and the hypergeomctric functions R±[ iv at the 
initial time of interest. (Note that because they oscillate as ln(Y) — > — oo, this initial 
time cannot be exactly zero.) To do this, one can go back to their defining differential 
equations, and solve them numerically, keeping only the solution that becomes a positive 
frequency plane wave as r — > +00. 

iii.c Superimpose on top of the classical Glasma fields (obtained in step ii), a linear combi- 
nation of the small fluctuations that are obtained by solving the eigenvalue equations in 
eqs. (97) and (104), multiplied by random Gaussian coefficients that have a variance 23 
given by eq. (65). 

iv. For each initial condition generated in this way, solve numerically the classical Yang-Mills 
equations forward in time, up to the proper time at which the observable should be evaluated. 
Repeat steps iii and iv in order to do a Monte-Carlo evaluation of the average over the 
fluctuations of the initial gauge fields. 



6 Applications and Outlook 

With this work, we have completed the resummation of all the leading contributions from quantum 
fluctuations to inclusive quantities at early times, in the Color Glass Condensate effective field 
theory approach to high energy nucleus-nucleus collisions. These quantum fluctuations can be 
factorized into v = and v ^ modes, where v is the Fourier conjugate to the space-time rapidity 
rj. The expression in eq. (5) described the contribution of v = modes, which correspond to 
summing the leading logarithmic a s ln(l/xi,2) contributions. Our final expression for the energy- 
momentum tensor, extending the expression in eq. (5) to include the leading secular terms, is 

(T^LLx+LInst. = J[D Pl Dp 2 }W Xl [pi}W X2 [p 2 ] 

x J[Va]F [a] T%[A[p u p 2 ] + a](x) , (106) 

where the weight functionals W Xli2 [pi,2] satisfy the JIMWLK renormalization group equation in 
eq. (6). The argument A = (A,£) denotes collectively the components of the classical fields and 
their canonicallly conjugate momenta on the initial proper time surface. These quantities are 
functionals of pi,p 2 and, as discussed previously, analytical expressions for these are available at 
t = + [47,51,55,56]. The initial spectrum of fluctuations Fo[a], defined in eq. (28), requires that 
one compute small fluctuations around the classical background fields A as r —¥ + . The formal 
expressions for these were derived in section 4. In section 5, we have outlined a practical algorithm 
to compute the path integral over small fluctuations. As numerical algorithms for solving the 
JIMWLK equation are now also available [65] and have been successfully implemented [66] , a full 
fledged numerical computation of eq. (106) is feasible in the near future. We should emphasize 
that this equation is valid for any inclusive quantity, not just the energy-momentum tensor. 

There are several applications of this formalism to understand key features of early time dy- 
namics in heavy ion collisions. We shall discuss a few of these here. 

23 Note that the (arbitrary) factor N K in this variance is cancelled by the normalization of the eigenfunctions a K 
and of the integration measure dfi K . Indeed, in eq. (64) one has c K ~ N K , a K ~ \/N K and dfi K ~ ^-/N K . 
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• Early thermalization It is important to emphasize that a numerical simulation of cq. (106) 
would describe the real time evolution of a quantum field theory that goes far beyond purely 
classical contributions but includes as well important quantum effects to all orders in per- 
turbation theory. (For similar considerations of the relative roles of classical versus quantum 
effects within the 2PI framework, we refer readers to refs. [67-69] and references therein.) In 
ref. [1], we developed the corresponding formalism for a scalar 4 theory, and demonstrated 
that the system developed an equation of state, allowing one to write the energy-momentum 
conservation equation for the resummed energy momentum tensor as the closed form set 
of equations corresponding to the ideal hydrodynamics satisfied by a perfect relativistic fluid. 
We interpreted this result as arising from a phase decoherence of the classical trajectories 
in eq. (106), with the different initial conditions given by the spectrum of initial quantum 
fluctuations. As discussed previously, we expect the same to occur in QCD as well. 

It is widely believed that the strong hydrodynamic behavior seen in heavy ion collisions 
requires early thermalization. However, as shown for the scalar theory, and may likely also 
be true for gauge theories, this is not a necessary condition for (nearly) perfect fluidity. 
It is interesting to ask what is the proper condition for thermalization in a quantum field 
theory. One such criterion is based on Berry's conjecture [40,41,43,44], that states that in a 
quantum system whose classical counterpart is chaotic, high lying energy eigenfunctions look 
like Gaussian random functions. It was later argued by Srcdnicki [42] that such eigenstates 
would appear to be thermal, for example, if one measures the single particle distribution. 
Since the underlying classical Yang-Mills theory is believed to be chaotic [45,46], the quantum 
system described by cq. (106) is a good candidate to explore these ideas in a quantum field 
theory. In our approach, the initial state is not an energy eigenstate, but rather a coherent 
state. It is formulated as a sum of plane waves with random Gaussian coefficients (see 
eqs. (64) and (65)) but there is no constraint that restrict the Gaussian random eigenstates 
to states of a given energy, in contrast to the states postulated by Berry. The self-interactions 
of the fields lead to a loss of this initial coherence, and it would be very interesting to see if 
thermalization occurs on the same time-scales as decoherence. 

Alternatively, one can look at the spectral function, defined in terms of the imaginary part 
of the retarded Green function, for the appearance of quasi-particle behavior, which allows 
for a kinetic theory description in terms of single particle distributions. It is interesting to 
explore whether there is a region of overlap between a description in terms of high occupation 
number fields and a kinetic theory description in terms of quasi-particles [70,71]. Such a 
regime may equivalently be described by a coupled set of equations for the classical fields and 
quasi-particle modes [72] which is reminiscent of the description of superfluids in condensed 
matter physics. Indeed, it is conceivable that the best description of such an overpopulated 
system, for intermediate times, might be as a Bose-Einstein superfluid [73], before the inelastic 
processes in the Glasma begin to dominate [74,75] and lead eventually to a conventional 
kinetic description. 

While the mechanism for thermalization is non-perturbative, it is still unclear whether the 
mechanism is weak coupling or strong coupling. One might anticipate that the former is more 
likely at higher energies due to the increasing dominance of semi-hard scales, but where such a 
transition occurs is unknown. It is intriguing that while the the mechanism of thermalization 
in AdS/CFT inspired strong coupling approaches appears very different (more "top down" 
than "bottom up"), there appear to be technical similarities between aspects of our weak 
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coupling approach and these approaches [76,77]. 

• Sphaleron transitions 

As emphasized, our master formula in cq. (106) is valid not only for the expectation value 
of the energy-momentum tensor but also for other inclusive quantities. One such example 
is the sphaleron transition rate r sp h a i. , which controls the mean squared change in the axial 
charge in thermal equilibrium through the relation 

((AQi,,) 2 )therm. = WtT sphal . , (107) 

where V is the spatial volume and q denotes quark flavor. In rcf. [78], it was noted that 
sphaleron transitions are not allowed in the boost invariant classical Glasma. Because the 
classical dynamics is effectively 2+1-dimensional, the second homotopy group of SU(3) gauge 
theory is zero, disallowing integer valued topological transitions. The quantum fluctuations 
we have been describing are no longer boost invariant, thereby allowing sphaleron transitions 
to go. In the non-equilibrium Glasma, the relevant quantity for computing the mean square 
change in the axial charge is the Wightman function [79] 

Glpi^Y) = ^F^Fr(X) J^F^FfiY)^ , (108) 

with F^ v = ^e^ a ^ F a p. The sphaleron transition rate is defined in terms of this quantity as 

r sph ai = J d A X G> P (X,0) . (109) 

This quantity may be computed by a formula similar to eq. (106). It has aroused much 
interest recently because in semi-peripheral heavy ion collisions, the combined effect of large 
external electromagnetic B fields and a large rate for topological transitions can lead to an in- 
duced charge separation phenomenon, called the Chiral Magnetic effect [80] , with observable 
consequences. Our approach allows in principle an ab initio computation of this effect. 

• Jet quenching 

Besides the large flow observed in heavy ion collisions, the apparent strong modification 
of rare final states (such as jets, high p± hadrons and heavy flavors) by the medium is 
adduced as confirmation of the high degree of opacity in the medium consistent with a 
strongly correlated fluid. The standard energy loss mechanism is radiative energy loss, which 
is dominated by collinear splittings that are primarily influenced by late time dynamics in a 
strongly interacting quark-gluon plasma. For a sampling of reviews, see refs. [81-84]. There 
are however potentially large modifications of the spectra of hard final states from energy 
loss at early times, that are are not included in this energy loss scenario [85]; for a recent 
treatment of large angle contributions in a framework similar to ours, see ref. [86]. However, 
no computation has thus far included next-to-leading order corrections to, for example, the 
single inclusive parton spectrum at early times taking into account both small x evolution 
and multiple scattering effects. This is done in our formalism when T^ v in cq. (106) is 
replaced by E p ^^] in fact, cq. (106) sums up a class of all order graphs that correspond 
to coherent multiple emissions where the momenta of all but one of the final state gluons 
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Figure 3: Example of graph included via the resummation of eq. (106). The green and red dots 
represent the color charges that describe the gluon content of the colliding nuclei in the CGC 
framework. 

is integrated over. A typical contribution is illustrated in the figure 3. Note that this 
diagram, corresponding to the resummed case, looks like a parton shower interacting with 
the background Glasma fields. However, unlike vacuum showers which can in space-time 
be visualized as being logarithmically divergent in the proper time, these "showers" are a 
consequence of the exponentially growing contributions from leading instabilities at each 
order in perturbation theory. 

We note that the single inclusive gluon spectrum also potentially suffers from collincar and 
infrared singularities. In the Glasma, it is possible that these could be regulated by strong 
multiple scattering and/or screening effects, but that remains to be proved. To avoid such 
complications, as in usual jet physics [87,88], one can look instead at correlators of the energy- 
momentum tensor that correspond to energy flow and are manifestly infrared and collinear 
safe. Clearly, there are a number of issues that need to be resolved here; the promising aspect 
of our formalism is that parton evolution, radiation and re-scattering can likely be treated, 
without ad hoc assumptions, in a consistent manner at early times. 
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A Alternative basis for the free field solution 

We will show in this appendix, that the form of the vacuum fluctuations can be written in terms 
of simple analytic functions if we trade the index v for an index 9 introduced via the following 
transformation 

dv 



ivO „H 



2tt 



a£(x) 



(110) 



where we denote k = (9,k±). 

To see that we are free to make this change in the quantum numbers, let us start from the 
real-valued field fluctuation as written in terms of the v coordinate 
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Then, making the trade from v to 9 we obtain 
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where we defined 



*kA 



dvcu_xe 



and in keeping with the notation employed throughout the text we have used k = (u, kj_) and 
k = (6, kj_). Since CkA is a Gaussian random variable so is its Fourier Transform, c?^ A , and we are 
free to generate random fluctuations using either basis. 

Using this new basis, the first physical solutions transforms to 
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After the transformation to the 9 variable, eqs. (82) can be integrated over time 
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(112) 



Note that eqs. (81) only specify a 2 ± and ot 2rj up to a r independent function. The r independent 
functions that can be added to a 2 ± and a 2v are not independent, because this modification must 
correspond to a residual gauge transformation. Thus the allowed modifications are 



a 2 j_ -> a 2 ±_ + f , a 2n -> a2>, + «(<9<?/) , 



(113) 
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where / is an arbitrary r independent function. We can choose to fix this residual gauge freedom 
by choosing the function / such that a^ 2 ( x ) vanishes at r = + . After fixing the residual gauge 
freedom, we have the following expression for the second physical solution 




a £ 2 ( x ) = __ | k,,!! | , "'■ (114) 

where we denote, 

<7_i_ = i tanh(0 — rj) 



2 ^ik±T cosh(6l — rj) 

1 _ e tk ± r cosh(e- v ) + ik ±TC osh(9 - rj) 
cosh 2 (6 — rj) 



01 = ~ • ( 115 ) 



B Wightman functions for free fields 

In this appendix, we shall construct the equal r Wightman function for free fields. The equal-time 
Wightman function corresponding to the first physical solution a£ is given by 

Of = J ^^^(r.^iJ^^Vyi) • (H6) 
Using the explicit expression in cq. (110), G^" can be written in the following formal way, 

G ? - <sp ( 4 -t ) * ( ia - 2 v c ^) 5 c - (Ajiibj) ■ ,117) 

where A is an infrared cutoff 24 and we have defined 

1 r rlr 

Ji(i/,d±) = - / d v e *^ e -tedxBmh(,/2) _ (llg) 

2 7 v 1 + a; 2 

Likewise, the Wightman function corresponding to the second physical solution is given by 

G *" = / f^K^'^)^'^)- ( 119 ) 

From the discussion following eq. (81), we noted that there was a residual gauge freedom remaining 
even after finding the two physical solutions. In the appendix A, we chose to fix this gauge by 
requiring that a v — > at r = + in order that Hamilton's equations are regular at r = 0. Clearly, 
for the problem at hand, this is the correct choice. However, computing the Wightman function 
will be made much easier by choosing a different gauge, 

I + ik ± r cosh(6 -r)) 

9± = * tanh(0 - rj) , g v = -5— . 

cosh (v — rj) 



The final result does not depend on this cutoff, thanks to the derivatives acting on the logarithm. 



33 



We should stress that there is a residual gauge freedom remaining in any correlator of gauge holds. 
With this new gauge choice we find 

^-(SfC 4 %" )^( ia " 2 \ /= ^) {< "-'' )b (Aklbii) < <120) 



with 



. . If dx 



2 cosh r) 



1 + 2a; 2 + cosh?? 



gii/rj e -ixd± sinh(r;/2) (121) 



We can also calculate the correlation function of the canonical momenta defined by e 1 = rd T ai. 
The free Wightman function for the transverse components of the electric held is defined as 

H iJ =J2 t2 J |S> [d T a£ x (T, V ,x ± )] [d T a^(T,rf,x' ± j\ . (122) 



A=l,2 



Note that the r derivatives remove the residual gauge freedom that still remained in the expression 
for G y as written above. This is expected since we don't expect any residual gauge freedom to 
remain when computing correlators of physical quantities such as th electric field. The Wightman 
function for the first physical solution is 

ffp =(^( it -i )*( i ^\f^) s ^-«'Hwi^r\) • (123) 

and similarly for the second physical solution 

ff? = -<^( o%, %' )^( i8 - 2 \/ I ^) { c-'') 1 °(Ajiibj) • (124) 

where the functions are defined as 

1 f dr 

F 3 (v, d±) = - - dr] [1 + 2x 2 + cosh??] e wr] e~ ixdl - sinh (»/2) ; ( 125 ) 

4 J VI + x 2 

1 f dr 

T A {v, d ± ) = -J - j= d V [1 + 2x 2 - cosh??] e w " e"^ sinh (»/2) . (126) 
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